# Подгружаем необходимые библиотеки import pandas as pdimport numpy as npfrom scipy import statsimport matplotlib.pyplot as pltimport itertoolsimport datetimefrom mpl_toolkits.basemap import Basemapimport foliumimport folium.plugins%matplotlib inline### Кратко воспроизводим итоги прошлой работы для получения необходимых данных# Загружаем файлdata = pd.read_csv('yellow_tripdata_2016-05.csv', sep=',', usecols = ['tpep_pickup_datetime', 'tpep_dropoff_datetime', 'passenger_count', 'trip_distance', 'pickup_longitude', 'pickup_latitude'])# Фильтруем по расстоянию:data = data[data.trip_distance>0]# Фильтруем по числу пассажировdata = data[data.passenger_count>0]# Фильтруем по длине поездкиdata = data[data.tpep_pickup_datetime < data.tpep_dropoff_datetime]# Фильтруем по принадлежности к городу lat_south = 40.49612lat_north = 40.91553long_west = -74.25559long_east = -73.70001data = data[(lat_south<data.pickup_latitude) & (data.pickup_latitude<lat_north) &(long_west<data.pickup_longitude)& (data.pickup_longitude<long_east)]# грузим сначала информацию по регионамregions = pd.read_csv('regions.csv', sep = ';')# Определяем интервалыbin_lat = np.unique(np.concatenate([np.unique(regions.south), np.unique(regions.north)]))bin_long = np.unique(np.concatenate([np.unique(regions.east), np.unique(regions.west)]))# Запускаем саму функциюregion_search = stats.binned_statistic_2d(data.pickup_latitude.values,data.pickup_longitude.values, None, 'count', bins=[bin_lat, bin_long], expand_binnumbers=True)# Сам регион можно получить из его групп по широте и долготе: data['region'] = region_search.binnumber[0]+50*region_search.binnumber[1]-50# Теперь добавим столбец со временем начала поездки с точностью до часаdata['hour'] = data.tpep_pickup_datetime.apply(lambda s: datetime.datetime.strptime(s[:13], '%Y-%m-%d %H'))# Далее нам потребуются только столбцы hour и regiondata = data[['hour','region']]# Группируем данные по часу и регионуgrouped_1 = data.groupby(['hour', 'region']).size()grouped_1.name = 'count'# Сформируем список всех часов в месяцеstart = datetime.datetime(2016, 5, 1)end = datetime.datetime(2016, 6, 1)hours = np.arange(start, end, datetime.timedelta(hours=1))# Формируем сетку час - регионproduct = itertools.product(hours, regions['region'].unique())grid = pd.DataFrame(list(product))grid.columns = ['hour', 'region']grid.set_index(['hour', 'region'], inplace=True)# Заполняем те позиции в сетке, которые отсутствуют в grouped_1 нулямиgrouped = grid.join(grouped_1)grouped.fillna(0, inplace=True)grouped['count'] = grouped['count'].astype(int)grouped.reset_index(inplace=True) # Теперь мы можем рассчитать количество регионов с нулевыми поездками, # а также количество регионов со средним числом поездок более 5trips_by_region = grouped.groupby('region')['count'].agg({'count': np.sum})print 'Количество регионов c нулевыми поездками: {}'.format((trips_by_region['count'] == 0).sum())# Рассчитаем среднее число поездок и количество регионов, для которых оно больше 5average_by_region = grouped.groupby('region')['count'].mean()print 'Количество регионов c числом поездок не менее 5: {}'.format(average_by_region[average_by_region>=5].count())# Собираем информацию по регионам для построения картregion_centres = []for r in regions.index: region = regions.loc[r, 'region'] mean = average_by_region[region] if mean >= 5: x = (regions.loc[r, 'west']+regions.loc[r, 'east'])/2 y = (regions.loc[r, 'south']+regions.loc[r, 'north'])/2 region_centres.append([y, x, mean])# Построим статическую карту с Empire State Buildingesb_longitude = -73.9857012esb_latitude = 40.7488006plt.figure(figsize=(15, 15))map = Basemap(llcrnrlon=long_west,llcrnrlat=lat_south,urcrnrlon=long_east,urcrnrlat=lat_north)map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 2000)map.plot(esb_longitude, esb_latitude, 'ok', markersize=12, color = 'r')plt.show()# Вспомогательные расчёты для статической картыx = np.linspace(long_west, long_east, 51)y = np.linspace(lat_south, lat_north, 51)x, y = np.meshgrid(x, y)array = np.array(trips_by_region['count'].values)z = np.reshape(array, (50, 50)).T# Сама статическая карта с иллюстрацией числа поездокplt.figure(figsize=(15, 15))static_map = Basemap(llcrnrlon=long_west,llcrnrlat=lat_south,urcrnrlon=long_east,urcrnrlat=lat_north)static_map.arcgisimage(service='ESRI_StreetMap_World_2D', xpixels = 2000)pcm = static_map.pcolormesh(x, y, z, alpha=0.4, cmap='RdBu')static_map.colorbar(pcm, 'bottom')plt.show()# Сформируем набор вспомогательных координатlib_lat = 40.6892494lib_long = -74.0466891map_centre = [(lat_south+lat_north)/float(2), (long_east+long_west)/float(2)]# Интерактивная карта со статуейinteractive_map = folium.Map(location=map_centre, zoom_start=12)folium.Marker([lib_lat, lib_long], popup='Statue of Liberty').add_to(interactive_map)interactive_map# Интерактивная карта с регионами с числом поездок по часамinteractive_map = folium.Map(location=map_centre, zoom_start=12)folium.plugins.HeatMap(region_centres, min_opacity=0.05, radius=5, blur=5, gradient={1: 'red'}).add_to(interactive_map)interactive_map